In seismic survey data, surface waves typically dominate intended reflection signals or body wave signals from the subsurface. Thus, it is desirable to attenuate them or remove them for further seismic processing. Current mitigation techniques typically assume the properties of the medium that transmits the surface waves are spatially homogeneous, often resulting in less than optimal surface wave mitigation and/or unwanted attenuation of reflection signals.
FIG. 1 depicts a typical process 100 to mitigate surface waves. Various existing filtering techniques may be used in the process 100. The process starts with one or more seismic records 101 for a particular region of interest. In block 102, the records 101 are analyzed at very few locations. The analysis involves determining the velocities and dispersion curves at the very few, selected locations. The data resulting from block 102 are sparsely sampled surface-wave properties 103. With this data 103, the process then designs filtering criteria to separate surface waves from body waves in block 104. The resulting sparse filtering criteria 105 are then interpolated by the process in block 106 for every location in the record and for every record in the input data 101. The interpolated criteria 107 are next used in block 108 by the process for the mitigation of surface waves in the input data 101 to produce data 109 with mitigated surface waves. Note that in other processes, the surface-wave properties are interpolated for every record instead of the filtering criteria being interpolated, but they result in the same inexact knowledge of the surface wave properties and/or the filtering criteria to separate surface waves from body waves.
One process that is used to reduce the effects of surface waves is phase-matched filtering, which is a method of removing the dispersion characteristics of the surface waves by flattening the surface waves in a seismic record. Phase matching also compresses the long and ringy surface-wave waveform in the time domain by removing the frequency-dependent velocity structure of the surface wave. This produces a surface wave that is not only flat but compact in the time-space domain of the seismic record. This compression of the surface wave is very advantageous because it allows small windows to be used over the limited frequency range of the surface wave to remove the surface wave. In an improvement to narrow time windows, Kim, U.S. Pat. No. 5,781,503, which is hereby incorporated herein by reference, teaches the use of a spatial low-pass filter on the time-aligned and compressed surface-wave data.
In phase-matched filtering, compression and alignment of surface waves are achieved by phase conjugating the surface waves G(f) in the frequency domain using the estimated phase velocity {circumflex over (v)}p(f). The surface waves are compressed in the time-domain after the phase conjugation, since the temporal elongation of the waveforms due to dispersion is negated. The phase conjugated waveforms are then aligned at t=to by a time shift implemented by a linear phase shift in the frequency domain, followed by the inverse Fourier transform. This can be mathematically expressed asĝc(t,{circumflex over (k)}r)=∫G(f)exp i[−{circumflex over (k)}rr−2 πf(t−to)]df,  (1)where {circumflex over (k)}r=2 πf/{circumflex over (v)}p(f), r=|r−rs|, rs and r are the locations of the source and the receiver, ĝc(t,{circumflex over (k)}r) is the waveform phase-conjugated by the phase term φ(r,f)={circumflex over (k)}rr and then time shifted to t=to.
Despite the value of aligning and compressing the surface waves, and the value of the subsequent spatial low-pass filtering, it is still necessary with phase-matched filtering of any kind to perform an analysis of the dispersion curves of the surface waves. These dispersion curves, or frequency-dependent phase velocities, are traditionally analyzed on some representative records from around the survey area. Seismic processors then typically apply one dispersion curve, {circumflex over (v)}p(f) to one group of traces, and another curve to another group of traces. In other words, the horizontal wavenumber {circumflex over (k)}r in Eq. (1) does not change within the group of traces, and thus spatial change of {circumflex over (k)}r within the group of traces cannot be accounted for.
The removal of phase to align a wavefield is practiced in several areas of geophysics. For example, the '503 patent to Kim applies phase removal to the alignment of surface waves and teaches the use of a single dispersion function to phase match all the traces in a seismic record under consideration. In standard seismic processing, a normal moveout (NMO) function is applied to prestack seismic data to align body wave reflections in a seismic record, see O. Yilmaz, Seismic Data Processing, Society of Exploration Geophysicists, 1987, which is hereby incorporated herein by reference. Again, only a single NMO function is applied to each trace in the record to achieve this alignment, though of course this single function results in a different time correction at each trace because it removes the effect of source-receiver offset distance. Use of a single function may be appropriate in common midpoint (CMP) processing when it is proper to ignore structuring and anisotropy, i.e., when the beds are essentially isotropic and horizontal, because the reflection represented on each trace in the CMP record is presumed by the sorting of the data to sample the same subsurface point.
When structural complexity is involved, NMO is no longer suitable, and prestack migration must be applied. In migration, phase matching or alignment of reflections in a record is accomplished by calculating the traveltime to the reflector as a wave propagates through a complex, spatially varying velocity overburden, see Yilmaz, ibid. The traveltime computation involves a path integration over the portions of the subsurface through which the wave travels on its path to the reflector and back to the surface. However, this path integral is a scalar integration assuming one number, e.g. velocity, for each cell, which is a volume element in three dimensional space, or voxel, that describes the subsurface along the wave path. Generalizations for anisotropy exist in which the traveltime is computed from a more general, vector velocity field that incorporates velocity as a function of direction. For anisotropy, the direction of the wave through the voxel is combined with the directional aspects of the velocity field to arrive at traveltime for the wave to the reflector and back to the surface. In all of these cases, a single traveltime is applied to each trace for each wavefield being phase corrected. In the simpler cases, the traveltime is derived from a single function, namely standard surface-wave phase matching and NMO. In the more complex cases, such as migration, the traveltime is derived from a different function for each trace, by path integration.
Another application, in many respects identical to seismic migration, is time-reversed focusing. In this application, acoustic wavefields received by a receiver array are time-reversed and then re-emitted into a medium in order to focus or image individual source points in the medium, for example see M. Fink and C. Prada, “Acoustic time-reversal mirrors,” Inverse Problems 17, R1-R38, 2001, which is hereby incorporated herein by reference. Since time reversal is equivalent to reversal of the sign of the phase in the frequency domain, this is equivalent to phase removal in a mathematical sense. In time-reversal, however, waves are physically retransmitted from a receiver array, and phase removal is achieved by the waves propagating through the medium. This is inherently different from the present invention where received wavefields are artificially back-propagated through the medium using knowledge of the spatially-varying velocity field of the medium. Although computational time-reversal techniques exist where physical retransmission of the waves is not required, see for example A. J. Berkhout, “Pushing the limits of seismic imaging, Part II: Integration of prestack migration, velocity estimation, and AVO analysis,” 1997, Geophysics, 954-969, which is hereby incorporated herein by reference, they are similar to true-amplitude migration. Furthermore, their purpose is directed to imaging.
For surface waves, path integration over the portions of the subsurface through which the wave travels is also known, see for example R. Snieder, “3-D linearized scattering of surface waves and a formalism for surface wave holography,” Geophys. J. R. astr. Soc. 84, 581-605, 1986, which is hereby incorporated herein by reference. However, the path integration is mostly used for the forward modeling of surface waves. Furthermore, these forward modeling formulations again incorporate amplitude terms at the source and the receiver locations, trying to account for both amplitude and phase of the surface waves. Even when it is used for phase-matching, for example see Stevens and McLaughlin, “Optimization of surface wave identification and measurement,” Pure appl. Geophys. 158, 1547-1582, 2001, which is hereby incorporated herein by reference, it was used to facilitate the detection and identification of weak surface wave events. Note that the goal is better localization of seismic sources in space.